library(plyr)
library(ggplot2)
library(sjmisc)
library(ordinal)

# Appendix D

a <- densityplot(study[study$condition == 2,"area"], 
                 xlab = "User-Defined Boundary Area - Study 1", bw = 1000, pch = 20)
b <- densityplot(study[study$condition == 1,"area"], 
                 xlab = "Location-Based Convex Hull Area - Study 1", bw = 1000, pch = 20)
c <- densityplot(replication[replication$condition == 2,"area"], 
                 xlab = "User-Defined Boundary Area - Study 2", bw = 500, pch = 20)
d <- densityplot(replication[replication$condition == 1,"area"], 
                 xlab = "Location-Based Convex Hull Area - Study 2", bw = 100, pch = 20)

print(a, split = c(1, 1, 2, 2), more = TRUE)
print(b, split = c(2, 1, 2, 2), more = TRUE)
print(c, split = c(1, 2, 2, 2), more = TRUE)
print(d, split = c(2, 2, 2, 2))

# Appendix E

###############################
##### Visual Outliers #########
###############################

outliers1 <- which(study$i == 222 | study$i == 374 | study$i == 749 | study$i == 392 | study$i == 680)
outliers2 <- which(replication$i == 113)

white.coefs <- rbind(std_beta(lm(prop.white ~ ch2white + pctwhite.zip + pctwhite.cty, study[-outliers1,]), include.ci = TRUE),  std_beta(lm(prop.white ~ ch2white + pctwhite.zip + pctwhite.cty, replication[-outliers2,]), include.ci = TRUE))
black.coefs <- rbind(std_beta(lm(prop.black ~ ch2black + pctblack.zip + pctblack.cty, study[-outliers1,]), include.ci = TRUE),  std_beta(lm(prop.black ~ ch2black + pctblack.zip + pctblack.cty, replication[-outliers2,]), include.ci = TRUE))
asian.coefs <- rbind(std_beta(lm(prop.asian ~ ch2asian + pctasian.zip + pctasian.cty, study[-outliers1,]), include.ci = TRUE),  std_beta(lm(prop.asian ~ ch2asian + pctasian.zip + pctasian.cty, replication[-outliers2,]), include.ci = TRUE))
hisp.coefs <- rbind(std_beta(lm(prop.hisp ~ ch2hisp + pcthisp.zip + pcthisp.cty, study[-outliers1,]), include.ci = TRUE),  std_beta(lm(prop.hisp ~ ch2hisp + pcthisp.zip + pcthisp.cty, replication[-outliers2,]), include.ci = TRUE))

college.coefs <- std_beta(lm(prop.college ~ ch2college + pctcollege.zip + pctcollege.cty, replication[-outliers2,]), include.ci = TRUE)
senior.coefs <- std_beta(lm(prop.seniors ~ ch2seniors + pctseniors.zip + pctseniors.cty, replication[-outliers2,]), include.ci = TRUE)

coefs <- rbind(black.coefs, hisp.coefs, asian.coefs, white.coefs, college.coefs, senior.coefs)

names <- factor(c(rep(c("Black", "Hispanic", "Asian", "White"), each = 6), rep("College-Educated", 3), rep("Senior Citizens", 3)), c("Black", "Hispanic", "Asian", "White", "College-Educated", "Senior Citizens"))

measures <- factor(rep(c("User-Defined Boundary", "ZCTA", "County"), length.out = 30), c("User-Defined Boundary", "County", "ZCTA"))

stn <- c(rep(c("July 2015", "July 2016"), each = 3, length.out = 24), rep("July 2016", 6))

ses <- (coefs[,1] - coefs[,2])/1.96

df <- data.frame(coefs, ses, names, measures, stn)

limits <- aes(ymax = ci.low, ymin= ci.hi)

p <- ggplot(df, aes(y=beta, x=measures, group=names, colour=names)) + geom_pointrange(limits, position=position_dodge(width=.5)) + theme_bw() + geom_hline(yintercept = 0) + labs(x = "Measures", y = "Estimated Regression Coefficients", color = "Demographics") + scale_color_grey(start = .8, end = .4)

p + facet_grid(stn ~ .) + theme(panel.grid.minor = element_blank(), panel.grid.major = element_blank())

white.coefs <- rbind(std_beta(lm(prop.white ~ pointswhite + combinedwhite + attachwhite + chwhite + pctwhite.zip + pctwhite.cty, study[-outliers1,]), include.ci = TRUE),  std_beta(lm(prop.white ~ naivepopwhite + combinedpopwhite + attachpopwhite + chpopwhite + pctwhite.zip + pctwhite.cty, replication[-outliers2,]), include.ci = TRUE))
black.coefs <- rbind(std_beta(lm(prop.black ~ pointsblack + combinedblack + attachblack + chblack + pctblack.zip + pctblack.cty, study[-outliers1,]), include.ci = TRUE),  std_beta(lm(prop.black ~ naivepopblack + combinedpopblack + attachpopblack + chpopblack + pctblack.zip + pctblack.cty, replication[-outliers2,]), include.ci = TRUE))
hisp.coefs <- rbind(std_beta(lm(prop.hisp ~ pointshisp + combinedhisp + attachhisp + chhisp + pcthisp.zip + pcthisp.cty, study[-outliers1,]), include.ci = TRUE),  std_beta(lm(prop.hisp ~ naivepophisp + combinedpophisp + attachpophisp + chpophisp + pcthisp.zip + pcthisp.cty, replication[-outliers2,]), include.ci = TRUE))
asian.coefs <- rbind(std_beta(lm(prop.asian ~ pointsasian + combinedasian + attachasian + chasian + pctasian.zip + pctasian.cty, study[-outliers1,]), include.ci = TRUE),  std_beta(lm(prop.asian ~ naivepopasian + combinedpopasian + attachpopasian + chpopasian + pctasian.zip + pctasian.cty, replication[-outliers2,]), include.ci = TRUE))

college.coefs <- std_beta(lm(prop.college ~ naivepctcollege + combinedpctcollege + attachpctcollege + chpctcollege + pctcollege.zip + pctcollege.cty, replication[-outliers2,], subset = naivepctcollege != Inf), include.ci = TRUE)

senior.coefs <- std_beta(lm(prop.seniors ~ naivex65plus + combinedx65plus + attachx65plus + chx65plus + pctseniors.zip + pctseniors.cty, replication[-outliers2,]), include.ci = TRUE)

coefs <- rbind(black.coefs, hisp.coefs, asian.coefs, white.coefs, college.coefs, senior.coefs)

names <- factor(c(rep(c("Black", "Hispanic", "Asian", "White"), each = 12), rep("College-Educated", 6), rep("Senior Citizens", 6)), c("Black", "Hispanic", "Asian", "White", "College-Educated", "Senior Citizens"))

measures <- factor(rep(c("Naive", "Time-Weighted", "Attachment-Weighted", "Convex Hull", "ZCTA", "County"), length.out = 60), c("Naive", "Time-Weighted", "Attachment-Weighted", "Convex Hull", "County", "ZCTA"))

stn <- c(rep(c("July 2015", "July 2016"), each = 6, length.out = 48), rep("July 2016", 12))

ses <- (coefs[,1] - coefs[,2])/1.96

df <- data.frame(coefs, ses, names, measures, stn)

limits <- aes(ymax = ci.low, ymin= ci.hi)

p <- ggplot(df, aes(y=beta, x=measures, group=names, colour=names)) + geom_pointrange(limits, position=position_dodge(width=.5)) + theme_bw() + geom_hline(yintercept = 0) + labs(x = "Measures", y = "Estimated Regression Coefficients", color = "Demographics") + scale_color_grey(start = .8, end = .4)

p + facet_grid(stn ~ .) + theme(panel.grid.minor = element_blank(), panel.grid.major = element_blank())

############################
##### 25th and 75th ########
############################

outliers1 <- which(study$area < quantile(study$area, .25, na.rm = T) | 
                     study$area > quantile(study$area, .75, na.rm = T)) 

outliers2 <- which(replication$area < quantile(replication$area, .25, na.rm = T) | 
                     replication$area > quantile(replication$area, .75, na.rm = T)) 

white.coefs <- rbind(std_beta(lm(prop.white ~ ch2white + pctwhite.zip + pctwhite.cty, study[-outliers1,]), include.ci = TRUE),  std_beta(lm(prop.white ~ ch2white + pctwhite.zip + pctwhite.cty, replication[-outliers2,]), include.ci = TRUE))
black.coefs <- rbind(std_beta(lm(prop.black ~ ch2black + pctblack.zip + pctblack.cty, study[-outliers1,]), include.ci = TRUE),  std_beta(lm(prop.black ~ ch2black + pctblack.zip + pctblack.cty, replication[-outliers2,]), include.ci = TRUE))
asian.coefs <- rbind(std_beta(lm(prop.asian ~ ch2asian + pctasian.zip + pctasian.cty, study[-outliers1,]), include.ci = TRUE),  std_beta(lm(prop.asian ~ ch2asian + pctasian.zip + pctasian.cty, replication[-outliers2,]), include.ci = TRUE))
hisp.coefs <- rbind(std_beta(lm(prop.hisp ~ ch2hisp + pcthisp.zip + pcthisp.cty, study[-outliers1,]), include.ci = TRUE),  std_beta(lm(prop.hisp ~ ch2hisp + pcthisp.zip + pcthisp.cty, replication[-outliers2,]), include.ci = TRUE))

college.coefs <- std_beta(lm(prop.college ~ ch2college + pctcollege.zip + pctcollege.cty, replication[-outliers2,]), include.ci = TRUE)
senior.coefs <- std_beta(lm(prop.seniors ~ ch2seniors + pctseniors.zip + pctseniors.cty, replication[-outliers2,]), include.ci = TRUE)

coefs <- rbind(black.coefs, hisp.coefs, asian.coefs, white.coefs, college.coefs, senior.coefs)

names <- factor(c(rep(c("Black", "Hispanic", "Asian", "White"), each = 6), rep("College-Educated", 3), rep("Senior Citizens", 3)), c("Black", "Hispanic", "Asian", "White", "College-Educated", "Senior Citizens"))

measures <- factor(rep(c("User-Defined Boundary", "ZCTA", "County"), length.out = 30), c("User-Defined Boundary", "County", "ZCTA"))

stn <- c(rep(c("July 2015", "July 2016"), each = 3, length.out = 24), rep("July 2016", 6))

ses <- (coefs[,1] - coefs[,2])/1.96

df <- data.frame(coefs, ses, names, measures, stn)

limits <- aes(ymax = ci.low, ymin= ci.hi)

p <- ggplot(df, aes(y=beta, x=measures, group=names, colour=names)) + geom_pointrange(limits, position=position_dodge(width=.5)) + theme_bw() + geom_hline(yintercept = 0) + labs(x = "Measures", y = "Estimated Regression Coefficients", color = "Demographics") + scale_color_grey(start = .8, end = .4)

p + facet_grid(stn ~ .) + theme(panel.grid.minor = element_blank(), panel.grid.major = element_blank())

white.coefs <- rbind(std_beta(lm(prop.white ~ pointswhite + combinedwhite + attachwhite + chwhite + pctwhite.zip + pctwhite.cty, study[-outliers1,]), include.ci = TRUE),  std_beta(lm(prop.white ~ naivepopwhite + combinedpopwhite + attachpopwhite + chpopwhite + pctwhite.zip + pctwhite.cty, replication[-outliers2,]), include.ci = TRUE))
black.coefs <- rbind(std_beta(lm(prop.black ~ pointsblack + combinedblack + attachblack + chblack + pctblack.zip + pctblack.cty, study[-outliers1,]), include.ci = TRUE),  std_beta(lm(prop.black ~ naivepopblack + combinedpopblack + attachpopblack + chpopblack + pctblack.zip + pctblack.cty, replication[-outliers2,]), include.ci = TRUE))
hisp.coefs <- rbind(std_beta(lm(prop.hisp ~ pointshisp + combinedhisp + attachhisp + chhisp + pcthisp.zip + pcthisp.cty, study[-outliers1,]), include.ci = TRUE),  std_beta(lm(prop.hisp ~ naivepophisp + combinedpophisp + attachpophisp + chpophisp + pcthisp.zip + pcthisp.cty, replication[-outliers2,]), include.ci = TRUE))
asian.coefs <- rbind(std_beta(lm(prop.asian ~ pointsasian + combinedasian + attachasian + chasian + pctasian.zip + pctasian.cty, study[-outliers1,]), include.ci = TRUE),  std_beta(lm(prop.asian ~ naivepopasian + combinedpopasian + attachpopasian + chpopasian + pctasian.zip + pctasian.cty, replication[-outliers2,]), include.ci = TRUE))

college.coefs <- std_beta(lm(prop.college ~ naivepctcollege + combinedpctcollege + attachpctcollege + chpctcollege + pctcollege.zip + pctcollege.cty, replication[-outliers2,], subset = naivepctcollege != Inf), include.ci = TRUE)

senior.coefs <- std_beta(lm(prop.seniors ~ naivex65plus + combinedx65plus + attachx65plus + chx65plus + pctseniors.zip + pctseniors.cty, replication[-outliers2,]), include.ci = TRUE)

coefs <- rbind(black.coefs, hisp.coefs, asian.coefs, white.coefs, college.coefs, senior.coefs)

names <- factor(c(rep(c("Black", "Hispanic", "Asian", "White"), each = 12), rep("College-Educated", 6), rep("Senior Citizens", 6)), c("Black", "Hispanic", "Asian", "White", "College-Educated", "Senior Citizens"))

measures <- factor(rep(c("Naive", "Time-Weighted", "Attachment-Weighted", "Convex Hull", "ZCTA", "County"), length.out = 60), c("Naive", "Time-Weighted", "Attachment-Weighted", "Convex Hull", "County", "ZCTA"))

stn <- c(rep(c("July 2015", "July 2016"), each = 6, length.out = 48), rep("July 2016", 12))

ses <- (coefs[,1] - coefs[,2])/1.96

df <- data.frame(coefs, ses, names, measures, stn)

limits <- aes(ymax = ci.low, ymin= ci.hi)

p <- ggplot(df, aes(y=beta, x=measures, group=names, colour=names)) + geom_pointrange(limits, position=position_dodge(width=.5)) + theme_bw() + geom_hline(yintercept = 0) + labs(x = "Measures", y = "Estimated Regression Coefficients", color = "Demographics") + scale_color_grey(start = .8, end = .4)

p + facet_grid(stn ~ .) + theme(panel.grid.minor = element_blank(), panel.grid.major = element_blank())


############################
###### 5th and 95th ########
############################

outliers1 <- which(study$area < quantile(study$area, .05, na.rm = T) | 
                     study$area > quantile(study$area, .95, na.rm = T)) 
outliers2 <- which(replication$area < quantile(replication$area, .05, na.rm = T) | 
                     replication$area > quantile(replication$area, .95, na.rm = T)) 


white.coefs <- rbind(std_beta(lm(prop.white ~ ch2white + pctwhite.zip + pctwhite.cty, study[-outliers1,]), include.ci = TRUE),  std_beta(lm(prop.white ~ ch2white + pctwhite.zip + pctwhite.cty, replication[-outliers2,]), include.ci = TRUE))
black.coefs <- rbind(std_beta(lm(prop.black ~ ch2black + pctblack.zip + pctblack.cty, study[-outliers1,]), include.ci = TRUE),  std_beta(lm(prop.black ~ ch2black + pctblack.zip + pctblack.cty, replication[-outliers2,]), include.ci = TRUE))
asian.coefs <- rbind(std_beta(lm(prop.asian ~ ch2asian + pctasian.zip + pctasian.cty, study[-outliers1,]), include.ci = TRUE),  std_beta(lm(prop.asian ~ ch2asian + pctasian.zip + pctasian.cty, replication[-outliers2,]), include.ci = TRUE))
hisp.coefs <- rbind(std_beta(lm(prop.hisp ~ ch2hisp + pcthisp.zip + pcthisp.cty, study[-outliers1,]), include.ci = TRUE),  std_beta(lm(prop.hisp ~ ch2hisp + pcthisp.zip + pcthisp.cty, replication[-outliers2,]), include.ci = TRUE))

college.coefs <- std_beta(lm(prop.college ~ ch2college + pctcollege.zip + pctcollege.cty, replication[-outliers2,]), include.ci = TRUE)
senior.coefs <- std_beta(lm(prop.seniors ~ ch2seniors + pctseniors.zip + pctseniors.cty, replication[-outliers2,]), include.ci = TRUE)

coefs <- rbind(black.coefs, hisp.coefs, asian.coefs, white.coefs, college.coefs, senior.coefs)

names <- factor(c(rep(c("Black", "Hispanic", "Asian", "White"), each = 6), rep("College-Educated", 3), rep("Senior Citizens", 3)), c("Black", "Hispanic", "Asian", "White", "College-Educated", "Senior Citizens"))

measures <- factor(rep(c("User-Defined Boundary", "ZCTA", "County"), length.out = 30), c("User-Defined Boundary", "County", "ZCTA"))

stn <- c(rep(c("July 2015", "July 2016"), each = 3, length.out = 24), rep("July 2016", 6))

ses <- (coefs[,1] - coefs[,2])/1.96

df <- data.frame(coefs, ses, names, measures, stn)

limits <- aes(ymax = ci.low, ymin= ci.hi)

p <- ggplot(df, aes(y=beta, x=measures, group=names, colour=names)) + geom_pointrange(limits, position=position_dodge(width=.5)) + theme_bw() + geom_hline(yintercept = 0) + labs(x = "Measures", y = "Estimated Regression Coefficients", color = "Demographics") + scale_color_grey(start = .8, end = .4)

p + facet_grid(stn ~ .) + theme(panel.grid.minor = element_blank(), panel.grid.major = element_blank())

white.coefs <- rbind(std_beta(lm(prop.white ~ pointswhite + combinedwhite + attachwhite + chwhite + pctwhite.zip + pctwhite.cty, study[-outliers1,]), include.ci = TRUE),  std_beta(lm(prop.white ~ naivepopwhite + combinedpopwhite + attachpopwhite + chpopwhite + pctwhite.zip + pctwhite.cty, replication[-outliers2,]), include.ci = TRUE))
black.coefs <- rbind(std_beta(lm(prop.black ~ pointsblack + combinedblack + attachblack + chblack + pctblack.zip + pctblack.cty, study[-outliers1,]), include.ci = TRUE),  std_beta(lm(prop.black ~ naivepopblack + combinedpopblack + attachpopblack + chpopblack + pctblack.zip + pctblack.cty, replication[-outliers2,]), include.ci = TRUE))
hisp.coefs <- rbind(std_beta(lm(prop.hisp ~ pointshisp + combinedhisp + attachhisp + chhisp + pcthisp.zip + pcthisp.cty, study[-outliers1,]), include.ci = TRUE),  std_beta(lm(prop.hisp ~ naivepophisp + combinedpophisp + attachpophisp + chpophisp + pcthisp.zip + pcthisp.cty, replication[-outliers2,]), include.ci = TRUE))
asian.coefs <- rbind(std_beta(lm(prop.asian ~ pointsasian + combinedasian + attachasian + chasian + pctasian.zip + pctasian.cty, study[-outliers1,]), include.ci = TRUE),  std_beta(lm(prop.asian ~ naivepopasian + combinedpopasian + attachpopasian + chpopasian + pctasian.zip + pctasian.cty, replication[-outliers2,]), include.ci = TRUE))

college.coefs <- std_beta(lm(prop.college ~ naivepctcollege + combinedpctcollege + attachpctcollege + chpctcollege + pctcollege.zip + pctcollege.cty, replication[-outliers2,], subset = naivepctcollege != Inf), include.ci = TRUE)

senior.coefs <- std_beta(lm(prop.seniors ~ naivex65plus + combinedx65plus + attachx65plus + chx65plus + pctseniors.zip + pctseniors.cty, replication[-outliers2,]), include.ci = TRUE)

coefs <- rbind(black.coefs, hisp.coefs, asian.coefs, white.coefs, college.coefs, senior.coefs)

names <- factor(c(rep(c("Black", "Hispanic", "Asian", "White"), each = 12), rep("College-Educated", 6), rep("Senior Citizens", 6)), c("Black", "Hispanic", "Asian", "White", "College-Educated", "Senior Citizens"))

measures <- factor(rep(c("Naive", "Time-Weighted", "Attachment-Weighted", "Convex Hull", "ZCTA", "County"), length.out = 60), c("Naive", "Time-Weighted", "Attachment-Weighted", "Convex Hull", "County", "ZCTA"))

stn <- c(rep(c("July 2015", "July 2016"), each = 6, length.out = 48), rep("July 2016", 12))

ses <- (coefs[,1] - coefs[,2])/1.96

df <- data.frame(coefs, ses, names, measures, stn)

limits <- aes(ymax = ci.low, ymin= ci.hi)

p <- ggplot(df, aes(y=beta, x=measures, group=names, colour=names)) + geom_pointrange(limits, position=position_dodge(width=.5)) + theme_bw() + geom_hline(yintercept = 0) + labs(x = "Measures", y = "Estimated Regression Coefficients", color = "Demographics") + scale_color_grey(start = .8, end = .4)

p + facet_grid(stn ~ .) + theme(panel.grid.minor = element_blank(), panel.grid.major = element_blank())

# Appendix F

cor(study[,c("ch2hisp", "pcthisp.cty", "pcthisp.zip")], use = "pairwise.complete.obs")[,1]
cor(study[,c("ch2white", "pctwhite.cty", "pctwhite.zip")], use = "pairwise.complete.obs")[,1]
cor(study[,c("ch2black", "pctblack.cty", "pctblack.zip")], use = "pairwise.complete.obs")[,1]
cor(study[,c("ch2asian", "pctasian.cty", "pctasian.zip")], use = "pairwise.complete.obs")[,1]

cor(replication[,c("ch2hisp", "pcthisp.cty", "pcthisp.zip")], use = "pairwise.complete.obs")[,1]
cor(replication[,c("ch2white", "pctwhite.cty", "pctwhite.zip")], use = "pairwise.complete.obs")[,1]
cor(replication[,c("ch2black", "pctblack.cty", "pctblack.zip")], use = "pairwise.complete.obs")[,1]
cor(replication[,c("ch2asian", "pctasian.cty", "pctasian.zip")], use = "pairwise.complete.obs")[,1]

cor(study[,c("combinedhisp", "pcthisp.cty", "pcthisp.zip")], use = "pairwise.complete.obs")[,1]
cor(study[,c("combinedwhite", "pctwhite.cty", "pctwhite.zip")], use = "pairwise.complete.obs")[,1]
cor(study[,c("combinedblack", "pctblack.cty", "pctblack.zip")], use = "pairwise.complete.obs")[,1]
cor(study[,c("combinedasian", "pctasian.cty", "pctasian.zip")], use = "pairwise.complete.obs")[,1]

cor(replication[,c("combinedpophisp", "pcthisp.cty", "pcthisp.zip")], use = "pairwise.complete.obs")[,1]
cor(replication[,c("combinedpopwhite", "pctwhite.cty", "pctwhite.zip")], use = "pairwise.complete.obs")[,1]
cor(replication[,c("combinedpopblack", "pctblack.cty", "pctblack.zip")], use = "pairwise.complete.obs")[,1]
cor(replication[,c("combinedpopasian", "pctasian.cty", "pctasian.zip")], use = "pairwise.complete.obs")[,1]

par(mfrow = c(2,2))

plot(full$pctwhite.zip, full$pointswhite, pch = 20, xlab = "% White (ZCTA)", ylab = "% White (Naive)", cex = .5)
lines(supsmu(full$pctwhite.zip, full$pointswhite), col = "red", lwd = 3)
abline(1, 1, lwd = 2)

plot(full$pctblack.zip, full$pointsblack, pch = 20, xlab = "% Black (ZCTA)", ylab = "% Black (Naive)", cex = .5, xlim = c(0,100), ylim = c(0, 100))
lines(supsmu(full$pctblack.zip, full$pointsblack), col = "red", lwd = 3)
abline(1,1)

plot(full$pctasian.zip, full$pointsasian, pch = 20, xlab = "% Asian (ZCTA)", ylab = "% Asian (Naive)", cex = .5, xlim = c(0,100), ylim = c(0, 100))
lines(supsmu(full$pctasian.zip, full$pointsasian), col = "red", lwd = 3)
abline(1,1)

plot(full$pcthisp.zip, full$pointshisp, pch = 20, xlab = "% Hispanic (ZCTA)", ylab = "% Hispanic (Naive)", cex = .5, xlim = c(0, 100), ylim = c(0, 100))
lines(supsmu(full$pcthisp.zip, full$pointshisp), col = "red", lwd = 3)
abline(1,1)

par(mfrow = c(2,2))

plot(full$pctwhite.zip, full$ch2white, pch = 20, xlab = "% White (ZCTA)", ylab = "% White (Boundary)", cex = .5)
lines(supsmu(full$pctwhite.zip, full$ch2white), col = "red", lwd = 3)
abline(1, 1, lwd = 2)

plot(full$pctblack.zip, full$ch2black, pch = 20, xlab = "% Black (ZCTA)", ylab = "% Black (Boundary)", cex = .5, xlim = c(0,100), ylim = c(0, 100))
lines(supsmu(full$pctblack.zip, full$ch2black), col = "red", lwd = 3)
abline(1,1)

plot(full$pctasian.zip, full$ch2asian, pch = 20, xlab = "% Asian (ZCTA)", ylab = "% Asian (Boundary)", cex = .5, xlim = c(0,100), ylim = c(0, 100))
lines(supsmu(full$pctasian.zip, full$ch2asian), col = "red", lwd = 3)
abline(1,1)

plot(full$pcthisp.zip, full$ch2hisp, pch = 20, xlab = "% Hispanic (ZCTA)", ylab = "% Hispanic (Boundary)", cex = .5, xlim = c(0, 100), ylim = c(0, 100))
lines(supsmu(full$pcthisp.zip, full$ch2hisp), col = "red", lwd = 3)
abline(1,1)

# Appendix G

white.coefs <- rbind(summary(clm(factor(seewhite) ~ ch2white + pctwhite.zip + pctwhite.cty, data = study))$coefficients[4:6,1:2], summary(clm(factor(seewhite) ~ ch2white + pctwhite.zip + pctwhite.cty, data = replication))$coefficients[4:6,1:2])
black.coefs <- rbind(summary(clm(factor(seeblack) ~ ch2black + pctblack.zip + pctblack.cty, data = study))$coefficients[4:6,1:2], summary(clm(factor(seeblack) ~ ch2black + pctblack.zip + pctblack.cty, data = replication))$coefficients[4:6,1:2])
hisp.coefs <- rbind(summary(clm(factor(seehisp) ~ ch2hisp + pcthisp.zip + pcthisp.cty, data = study))$coefficients[4:6,1:2], summary(clm(factor(seehisp) ~ ch2hisp + pcthisp.zip + pcthisp.cty, data = replication))$coefficients[4:6,1:2])
asian.coefs <- rbind(summary(clm(factor(seeasian) ~ ch2asian + pctasian.zip + pctasian.cty, data = study))$coefficients[4:6,1:2], summary(clm(factor(seeasian) ~ ch2asian + pctasian.zip + pctasian.cty, data = replication))$coefficients[4:6,1:2])

coefs <- rbind(black.coefs, hisp.coefs, asian.coefs, white.coefs)

names <- factor(c(rep(c("Black", "Hispanic", "Asian", "White"), each = 6)))

measures <- factor(rep(c("User-Defined Boundary", "ZCTA", "County"), length.out = 24))

stn <- c(rep(c("July 2015", "July 2016"), each = 3, length.out = 24))

ses <- coefs[,2]
coefs <- coefs[,1]

df <- data.frame(coefs, ses, names, measures, stn)

df$measures <- relevel(df$measures, "User-Defined Boundary")

limits <- aes(ymax = coefs + ses*1.96, ymin = coefs - ses*1.96)

p <- ggplot(df, aes(y=coefs, x=measures, group=names, colour=names)) + geom_pointrange(limits, position=position_dodge(width=.5)) + theme_bw() + geom_hline(yintercept = 0) + labs(x = "Measures", y = "Estimated Ordered Logit Coefficients", color = "Demographics", title = "User-Defined Boundary Models") + scale_color_grey(start = .8, end = .4)

p + facet_grid(stn ~ .) + theme(panel.grid.minor = element_blank(), panel.grid.major = element_blank())

white.coefs <- rbind(summary(clm(factor(seewhite) ~ pointswhite + combinedwhite + attachwhite + chwhite + pctwhite.zip + pctwhite.cty, data = study))$coefficients[4:9,1:2], summary(clm(factor(seewhite) ~ naivepopwhite + combinedpopwhite + attachpopwhite + chpopwhite + pctwhite.zip + pctwhite.cty, data = replication))$coefficients[3:8,1:2])

black.coefs <- rbind(summary(clm(factor(seeblack) ~ pointsblack + combinedblack + attachblack + chblack + pctblack.zip + pctblack.cty, data = study))$coefficients[4:9,1:2], summary(clm(factor(seeblack) ~ naivepopblack + combinedpopblack + attachpopblack + chpopblack + pctblack.zip + pctblack.cty, data = replication))$coefficients[4:9,1:2])

hisp.coefs <- rbind(summary(clm(factor(seehisp) ~ pointshisp + combinedhisp + attachhisp + chhisp + pcthisp.zip + pcthisp.cty, data = study))$coefficients[4:9,1:2], summary(clm(factor(seehisp) ~ naivepophisp + combinedpophisp + attachpophisp + chpophisp + pcthisp.zip + pcthisp.cty, data = replication))$coefficients[4:9,1:2])

asian.coefs <- rbind(summary(clm(factor(seeasian) ~ pointsasian + combinedasian + attachasian + chasian + pctasian.zip + pctasian.cty, data = study))$coefficients[4:9,1:2], summary(clm(factor(seeasian) ~ naivepopasian + combinedpopasian + attachpopasian + chpopasian + pctasian.zip + pctasian.cty, data = replication))$coefficients[4:9,1:2])

coefs <- rbind(black.coefs, hisp.coefs, asian.coefs, white.coefs)

names <- factor(c(rep(c("Black", "Hispanic", "Asian", "White"), each = 12)))

measures <- factor(rep(c("Naive", "Time-Weighted", "Attachment-Weighted", "Convex Hull", "ZCTA", "County"), length.out = 48))


stn <- c(rep(c("July 2015", "July 2016"), each = 6, length.out = 48))

ses <- coefs[,2]
coefs <- coefs[,1]

df <- data.frame(coefs, ses, names, measures, stn)

df$measures <- factor(df$measures, levels = c("Naive", "Time-Weighted", "Attachment-Weighted", "Convex Hull", "County", "ZCTA"))

limits <- aes(ymax = coefs + ses*1.96, ymin = coefs - ses*1.96)

p <- ggplot(df, aes(y=coefs, x=measures, group=names, colour=names)) + geom_pointrange(limits, position=position_dodge(width=.5)) + theme_bw() + geom_hline(yintercept = 0) + labs(x = "Measures", y = "Estimated Ordered Logit Coefficients", color = "Demographics", title = "Location-Based Estimate Models") + scale_color_grey(start = .8, end = .4)

p + facet_grid(stn ~ .) + theme(panel.grid.minor = element_blank(), panel.grid.major = element_blank())


